


####################################################
####          Replication Code for              ####
####     Spatial Interdependence                ####
####     and Instrumental Variable Models       ####
####################################################
####                                            ####
####                9/ 06/2018                  ####
####  this file runs R-file to create W's       ####
####      for Ramsay Application                ####
####                                            ####
###                                             ####
####################################################


data = read_dta("replicationOIL.dta")
head(data)
### table 3 model 4
### xi: ivreg2 normpolity2 (logoilrevpc=logORDdef2k) loggdppc gdpgrowth polity68 i.year, robust bw(2) first small rf endog(logoilrevpc)


### subset to complete cases
data_complete= data[complete.cases(data[,c("normpolity2", "logoilrevpc", "logORDdef2k", "loggdppc", "gdpgrowth", "polity68")]), ]
names(data_complete)


## now we need to add spatial data for each year 
## first add country codes

data_complete$COWCODE = countrycode(data_complete$country_name, "country.name", "cown", warn = TRUE)
data_complete = data_complete[order(data_complete$year,data_complete$COWCODE),]

years = unique(data_complete$year)


#### create dataframe that has the only countries for which we can get the shapefile
shapefile_countries = data.frame(data_complete[0,])
  
for(i in 1:length(years)){
  cshp = cshp(date=as.Date(paste(years[i],"-6-30", sep = "")), useGW=TRUE)
  dat = data_complete[data_complete$year == years[i],]
  dat = dat[dat$COWCODE %in% cshp$COWCODE,]
  shapefile_countries = bind_rows(shapefile_countries,dat)
}


#### function to create W matrix for each year, since countries vary by year
W_creator = function(year){
  cshp = cshp(date=as.Date(paste(year,"-6-30", sep = "")), useGW=TRUE)
  dat = data_complete[data_complete$year == year,]
  dat = dat[dat$COWCODE %in% cshp$COWCODE,]
  cshp = cshp[cshp$COWCODE %in% dat$COWCODE,]
  cshp@data = join(cshp@data, dat, by = "COWCODE", ) 
  cshp = cshp[order(cshp$COWCODE),] 

  neighbors = poly2nb(cshp)
  lw = nb2listw(neighbors, style="W", zero.policy=TRUE)
  W = listw2mat(lw)
  return(W)
}


#### create empty list, then fill with W matrix for each year
WW = list(NULL)
WW = lapply(1:length(years), function(x) W_creator(years[x]))
            
#### create large W matrix
WWW = bdiag(WW)
WWW = as(WWW,"matrix")
write.table(WWW, "Ramsay_W.csv", sep=",", col.names=F, row.names=F)

#### order data so that it fits the W matrix
shapefile_countries <- shapefile_countries[order(shapefile_countries$year, shapefile_countries$COWCODE),]
#### create ID variable for Stata
shapefile_countries$id = seq(1:dim(shapefile_countries)[1])
write_dta(shapefile_countries, "Ramsay_data_Table2.dta")
